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We investigate the nature of the parallel-roll to spiral-defect-chaos (SDC) transition in Rayleigh- 
Benard convection, based on the generalized Swift-Hohenberg model. We carry out extensive, 
systematic numerical studies by, on one branch, increasing the control parameter gradually from 
the parallel roll regime to the SDC regime and, on the other branch, decreasing it in the opposite 
manner. We find that the data of several time-averaged global quantities all form hysteretic loops 
from the two branches. We also discuss several possible scenarios for the transition and analyze our 
data for SDC accordingly. We conclude that the roll-to-SDC transition is first-order in character 
and that the correlation length diverges at the conduction to convection onset. We further postulate 
that this transition can be understood somewhat similar to the hexagon-to-roll transition in non- 
Boussinesq fluids. Finally we comment on the differences between our conclusion and those in two 
experiments. 
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I. INTRODUCTION 

Many nonequilibrium systems exhibit self-organized 
pattern- forming phenomena jlj. In the past few years, 
a new type of intrinsic pattern has emerged in many 
disciplines of science. These patterns are characterized 
by their extensive, irregular behavior in both space and 
time, which are known as spatiotemporal chaos (STC) 
0. STC typically exists in large size systems and its 
complexity increases dramatically with the system size 
[gj. Owing to its generic dynamical complexity, STC 
hence poses a great challenge to both experimentalists 
and theoreticians. 

From the very beginning, Rayleigh-Benard convection 
(RBC) has been a paradigm in the study of pattern for- 
mation in driven dissipative systems, because of its rel- 
ative simplicity and high precision in controlled experi- 
ments §]. RBC can occur when a thin horizontal layer 
of fluid is heated from below. The system is described 
by three dimensionless parameters (!]]: (a) the Rayleigh 
number R = gad 3 AT/ nv, in which g is the gravitational 
acceleration, d the layer thickness, AT the temperature 
gradient across the layer, a the thermal expansion co- 
efficient, k the thermal diffusivity and v the kinematic 
viscosity; (b) the Prandtl number a = v j 'k; and, (c) the 
aspect ratio T = L/2d where L is the horizontal size of 
the system. The Rayleigh number R is the control pa- 
rameter of the system; the Prandtl number a specifies 
the fluid properties. It is convenient to introduce a re- 
duced control parameter e = (R — R c )/R c , where R c is 
the critical value of R at which the fluid bifurcates from 
a static conductive state to a convective state. 

RBC has been studied extensively in the literature 
jl],|| ■ Theoretical analyses by Busse and his coauthors ||] 
predict that parallel roll states are stable inside a stability 
domain in (R, k, a) space with k the wavenumber, which 
is known as the "Busse balloon" . Surprisingly, recent ex- 



periments (HJ6) and numerical studies (7[|| , using systems 
with a ~ O(l) and large T, revealed that the parallel roll 
state yields to a spatiotemporally chaotic state even for 
states inside the Busse balloon. This spatiotemporally 
chaotic state, called spiral-defect-chaos (SDC), exhibits 
very complicated dynamics, illustrated by the interplay 
of numerous rotating spirals, patches of moving rolls, in- 
tricate grain boundaries, dislocations and other defects 
J^raj . Its discovery has since stimulated many experimcn- 
talt||Jl^-|l|, theoretical [[|Jl4| and numerical [0-f|l| 
efforts to understand it. 

The nature of the parallel roll to SDC transition is 
one of the important questions with respect to SDC and 
has been investigated in several experimental and nu- 
merical studies [ pl]Jl^Jl5| ] . By solving the generalized 
Swift-Hohenberg (GSH) model of RBC [||-f§ for non- 
Boussinesq fluids with random initial conditions, two of 
us []l5|| characterized the transition by the behavior of 
time-averaged global quantities such as convective cur- 
rent J, vorticity current Q (called vortex energy in Ref. 
p"5[) and spectra entropy 5 |l9[ ]. It was found that the 
convective current seems to be smooth across the transi- 
tion temperature ex but both the vorticity current and 
the spectra entropy seem to obey power-law behavior 
near ct- However this study was unable to distinguish 
between a gradual or sharp transition. Despite its incon- 
clusiveness on the nature of the roll-to-SDC transition, 
this study suggests that studying such time-averaged 
global quantities is quite useful. It hence motivated us 
to develop a phenomenological theory for STC, including 
SDC, in RBC @. In the theory, we made a random 
phase approximation for the spatiotemporally chaotic 
states and assumed that their time-averaged structure 
factor S(k) satisfies a scaling form with respect to the 
two-point correlation length £2- With these assumptions, 
we obtained analytical expressions for both J and Q in 
terms of measurable quantities. These theoretical results 
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provide us some new insights on the nature of SDC. In ad- 
dition, a recent experimental study also found spec- 
tra entropy 5 to be a useful quantity. 

On the experimental side, Morris et al. studied the 
structure of SDC using a circular cell. They found that 
the correlation length £2 is smooth across er and diverges 
at e = with a mean- field exponent. But the data for 
the correlation time r was consistent with either a diver- 
gence at £t with a mean-field exponent or a divergence at 
e = with a non-mean-field exponent. However, a differ- 
ent conclusion was reached by Cakmur et al. p2[ recently, 
who used a square cell in their experiment. They found 
that £2 diverges at ct with a small exponent. In critical 
phenomena, we know that finite (infinite) £2 and t are 
normally associated with first-order (second-order) tran- 
sitions. If this is also true in non-equilibrium phenomena, 
then these two experiments should lead to different con- 
clusions about the nature of the roll-to-SDC transition at 
ey. Although it is possible that the nature of the tran- 
sition is different for systems with different geometry, as 
suggested in Ref. |l2|], it is not clear whether this is true 
in the limit of an infinite system. We believe, however, 
that the results of Cakmur et al. |12| should be treated 
with caution, particularly given the absence of data for 
£2 for the parallel roll states (i.e., e < er) and the rela- 
tively limited number of data near er. We will comment 
further about these two experiments in Sec. 

In this paper, we present our extensive, systematic nu- 
merical studies of SDC, based on the GSH model of RBC 
for Boussinesq fluids. In comparison with studies in Ref. 
p5[ , we use random initial conditions at e = 0.05 and 
e = 0.8 only. Then, after completing the calculation at 
one e, we increase e (originally from e = 0.05) or decrease 
e (originally from e = 0.8) gradually and take the final 
state from the previous e as our initial condition. We 
hence obtain two different branches of data, one from 
increasing e and the other from decreasing e. We find 
that the results for £2, J, fi and H all form hysteretic 
loops from the two branches during the roll-to-SDC tran- 
sition. We analyze our data in accordance with our the- 
oretical results jL4|. We conclude that the roll-to-SDC 
transition is first-order in character. We also postulate 
that this transition can be understood somewhat similar 
to the hexagon-to-roll transition in non-Boussinesq flu- 
ids ^0|; namely, we postulate that the SDC bifurcation 
actually occurs at e = 0. But, since SDC is unstable 
(or metastable) against parallel roll states at smaller e, 
it emerges only for e > ex- 

This paper is organized as follows. In Sec. |l[ we 
introduce the GSH model of RBC and define some im- 
portant time-averaged global quantities. We then sum- 
marize our theoretical results in Ref. |14|. In Sec. [II, we 



ture of the roll-to-SDC transition and comment on the 
different conclusions between this work and the earlier 
experimental studies [ pl|^2[ . 

II. THEORETICAL RESULTS 



The GSH model of RBC |J-[18| is widely accepted 
for theoretical study. This model is derived from the 
three-dimensional hydrodynamic equations, but is much 
simpler to study both numerically and analytically. The 
GSH model contains two coupled equations in two- 
dimensional space r = (x, y), one for the order parameter 
ip(r,t) and the other for the mean-flow field ((r,t). The 
convective patterns in RBC are completely determined 
by the order parameter tp(r, t). The amplitude equations 
for the GSH model and the hydrodynamical equations 
are the same in the leading order near onset. Numerical 
solutions of this model or its modified versions have not 
only reproduced most patterns observed in experiments 
but also resembled experimental results relatively well 
]l|,[7U^ jr^,p^|,pT[ | . But there are some shortcomings in the 
model ]9|,p2|: The stability boundary of the model does 
not coincide with that of hydrodynamics; it induces an 
unphysical, short-ranged cross roll instability; and both 
the shape and the peak position of the power spectrum 
for SDC are different from those in the real system. Even 
so, owing to its simplicity and its qualitative resemblance 
to real systems, this model is very valuable in studying 
RBC. 

In the GSH model, the order parameter i/j(r, t) satisfies 



d t rjj + g m U • W = [e - (V 2 + l) 2 ] V - ^ 



(1) 



where V is the gradient operator in two-dimensions and 
U(r) is the mean-flow velocity given by U(r) = V£(r, t) x 
e z . The mean- flow field £(r, t), on the other hand, satis- 
fies g§ 

[d t - <r(V 2 - c 2 )] V 2 C - e 2 ■ [V(V 2 V) x W>] . (2) 

Variables in these equations have been rescaled for nu- 
merical convenience. Their relations to their physical 
values can be found in Ref. fli) . [See Eqs. (1) - (3) and 
(58) - (60) there. 1 For example, the reduced Rayleigh 
number e in Eq. ([[]) is related to its physical value e oxp t 
by e CX pt = 0.3594e. The rescaling factors for ip, £, t and 
r and parameters g m , a and c 2 can also be found in Ref. 

& 

We now define several important time-averaged global 
quantities in RBC. The first one is the time-averaged 
convective current defined as 

j = a- 1 f dr^M) = ^v;*(k,^(k,i), ( 3 ) 

where ^>(k, t) is the Fourier component of ^(r, t), A is the 
we discuss the subtleties involved in determining the na- area of the system and F(t) represents the time-average 



discuss possible scenarios with regard to SDC, according 
to whether the correlation length £2 diverges at e = or 
at ct- We present the details of our numerical studies 
in Sec. IV. We also analyze the data for £2, J, and 
S to test these different scenarios. In the last section, 
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of F(t). This quantity increases from J = to J > at 
the conduction-to-convection onset and hence character- 
izes the transition. The second one is the time-averaged 
vorticity current defined as 



n = A- 



druj 2 (r, t), 



(4) 



where u) z (r,t) — — V 2 £(r, t). This quantity reflects the 
distortion of patterns at large distance. It is identically 
zero for perfect parallel rolls and increases dramatically 
for SDC |l5| ]. Because of this, it has been speculated that 
one may take f2 as one of the order parameters to char- 
acterize the roll-to-SDC transition The third one is 
the time-averaged spectra entropy defined as [Tg| ] 

H=-^5(k,t)ln5(k,t), (5) 

k 

where the structure factor S*(k, t) is given by 

S(k,*)=^*(k,*)^(k,t)/J. (6) 

Like its counterpart in thermodynamics, the spectra en- 
tropy S is related to the randomness of all excited states. 
Its value is In 2 = 0.6931 for perfect parallel roll states 
but increases dramatically for SDC. This quantity was 
also found useful to characterize the roll-to-SDC transi- 
tion JlJJl^] • Another important quantity is the two-point 
correlation length defined as 



t2=[(k 2 ) k -(k)l] 



-1/2 



(7) 



where we has used the notation (F(k))k = 
J dk kS(k)F(k) in which S(k) is the time-averaged 
and azimuthally averaged structure factor normalized by 
J °° dkkS(k) = 1. Clearly £2 specifies the width of S(k). 

In Ref. |14| we presented our analytical calculations, 
using the GSH model, of J and Q for STC in RBC. These 
calculations are valid for both SDC and phase turbulence 
(PT) |2^,Q. By assuming the time-averaged two-point 
correlation function 



C(r l5 r 2 ) = ^(ri.tWfa.ty^ri,*), 



(8) 



is translation invariant in STC, i.e., C(ri,r2) = C(rj — 
r2), we found that the phases of two i/>(k, t) fields are mi- 
correlated in time unless they have the same wavenumber 
k. Furthermore, we applied a random phase approxima- 
tion (RPA) to STC in which four-point correlation func- 
tions are approximated by products of two-point correla- 
tion functions such as ipipipip ~ ipip ipip. Using this RPA, 
we derived J and Q in terms of S(k). We further as- 
sumed that the structure factor satisfies a scaling form 
with respect to £2, i-e., 



kS(k) = &T[(k 
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(9) 



where k max is the peak position of kS(k) and F{x) is the 
scaling function satisfying dx T(x) = 1. [Since k > 



in kS(k), the lower limit for 3-{x) is —k max £,2i which we 
approximate by —00.] From these assumptions, we ob- 
tained explicit formulas for both J and ft in the leading 
order of £j" . For SDC, these results can be written as 
El 



and 
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Hx 2 ) x 



e 2 



\' J SDC> 



(10) 



(11) 



where we have used the notation (F(x)) x — 
J™ dx T(x)F(x). Inserting k = k max + and Eq. 

(g) into Eq. (0), it is easy to see that {x 2 ) x = 1 + {x) 2 x > 1 
and (fc)fe = k m ax + €2 (x) x - 

In comparison, the convective current for perfect par- 
allel rolls with wavenumber fcg has been evaluated to be 
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J Roll 



§[e - (1 - fc ) 2 ] |17|. To our knowledge, there 

SH 



is no explicit formula for J%^ u for distorted rolls. If one 
uses this expression for but replaces (1 — fco) 2 with 
((1 — k) 2 )k to account for the finite width of the power 
spectrum, one finds for distorted roll states that 
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(12) 



where £2 and {x 2 ) x have exactly the same meanings as 
in Eq. ([l7j|), but their values may be different. 

Notice that the formulas for J§p C and JnSi are the 
same. This, however, is due to the simplification that 
the coupling constant of the nonlinear term, ip 3 , is taken 
as a constant in Eq. (|l|). In a more realistic description 
of hydrodynamics, this coupling constant, say qfcos a), 
is angle dependent and has been evaluated in Ref. |17f] 
[before the rescalings leading to Eq. (|J)]. For such a 
coupling constant, the time-averaged convective current 
for SDC, before the rescalings, has been calculated in 
Ref. [fli as 



JSDC 



9SDC 



(13) 



where £ 2 ~ 0.148 @ and g SDC = 1.1319 + 0.0483CT- 1 + 
0.0710er~ 2 jL4j. Correspondingly, the convective current 
for parallel rolls is 
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f 2 



with gjiou = 0.6995 - 0.0047CT- 1 + 0.0083cr~ 2 
these expressions, one finds that 



(14) 
From 



A J = Jsdc — JroI, 



1 



9sdc 



gSDC 9Roll 



SDC 



9Roll 



(15) 
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Since 2/gsDC is not equal to 1/gRoii for most values 
of cr, the first term above is not zero, ft is highly 
unlikely that this term may be cancelled by the con- 
tributions from the two other terms. Thus, there ex- 
ist discontinuities in the value and the slope of J dur- 
ing the roll-to-SDC transition. This is not surprising 
considering that J depends sensitively on the structure 
of the convective pattern |l7j and that the structures 
of parallel rolls and SDC are so different. Assuming 
[( x2 )x/^]sdc = [{x 2 ) x /(,1]roII at the transition temper- 
ature, we find that AJ/J Roll = 0.1239 for a = 1. So 
the value and the slope of J jump about 10% during the 
roll-to-SDC transition. Similar jumps can be found for 
other values of a under the same assumption; see Fig. 1. 

A.J °-3 



FIG. 1. I J Ron vs. 
AJ = Jsdc — JroII- 
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Prandtl number a where 



III. POSSIBLE SCENARIOS 

From Eqs. ( [To| ) and (|TT|) , it is obvious that the behav- 
ior of Jsdc and flsDC depend sensitively on the two- 
point correlation length £2- For simplicity, we drop the 
superscripts in Eqs. ( |l0|) - (|l2|). We assume that £2 has 
a power-law behavior such as 



£2 ~ £2,0^ 



(16) 



where e is the basic scaling field for SDC. Similar behav- 
ior has been found for £2 in phase turbulence (PT) f p"4[ . 
In that case, since the transition to PT occurs at e = 
P^pi} , one has simply e — e 113]. The possible scenario 
for SDC, however, is more subtle since the roll-to-SDC 
transition occurs at a positive temperature ex [||~@] ■ One 
obviously has two alternative choices for the scaling field 
einEq. ©: (A) e = e or (B)e = e-e T . The case (A) is 
similar to the situation in the hexagon-to-roll transition 
pO| where the two-point correlation length £2 is finite 
at the transition temperature ex but diverges at e = 0. 
This scenario is consistent with the experimental result 
by Morris et al. The case (B) resembles the situ- 

ation in critical phenomena in which £2 diverges at ex- 



This scenario was suggested by Cakmur et al. J12| . We 
now discuss these two scenarios separately. 

(A) If e = e. This implies that all properties of SDC 
are controlled at e = rather than at e — ex- As far as 
scaling relations are concerned, this case is similar to the 
situation in PT fl4|| . Similar to those in PT, one may 
hence define power laws such as 



Js 



DC 



and Qsdc ~ Sloe 



One finds from Eqs. 
relation 



(17) 



ll|) and (hq) the following scaling 



A = 2/i + 2v. 



(18) 



In comparison, one has A = 2/j, + v in PT |l4[ . From Eq. 
(]l0|), since Jsdc is positive by definition, the values of 
the exponents satisfy 



M= 1 , v>\/2 and A = 2 + 2v > 3. 



(19) 



Now it is useful to further distinguish two different 
cases: (Al) v = 1/2, or (A2) v > 1/2. Case (Al) cor- 
responds to a mean-field exponent v, in which one has 
that 



and 
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Jn = x 
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2^2 



1 - 



4(x 2 



f 2 

S2,0 



V4c 2 + ? 



- 1 



J 2 

£2,0 



(20) 



(21) 



The amplitudes Jo and Qq depend on two phenomeno- 
logical parameters £2,0 and (x 2 ) x . In case (A2), the ex- 
ponent v has a non-mean-field value. Now, since 2v > 1, 
the ^2~ 2 ~ e2v term in Eq. jl0| ) only adds a correction 
to the leading singularity. Instead of Eq. ( p"7j ) , one may 
define 



Jsdc = Joe^l + he" 



and 



n 



SDC 



= fV A [l + wie 



Ai 



(22) 



(23) 



with fix = 2v — fi > and Ai > 0. Consequently, instead 
of Eq. (20), one has that 



Jo = 2/3 and = -4(x 2 ) x /£ 



(24) 



While Slo is still given by Eq. ( |2l[ ), one must use the 
corresponding new value of Jo- The values of uj\ and 
Ai, however, cannot be determined without knowing the 
behavior of £2 beyond the leading term described in Eq. 

(B) If e = e — ex- As we mentioned, this case resembles 
the situation in critical phenomena in which £2 diverges 
at ex- Then, as in critical phenomena, other quantities 
should also have singular behaviors at ex- Now, instead 
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of Eq. (g^) or Eqs. and (J23J), we define, near I = 0+, 
that 



Jsdc ~ Joe - ^s,oe^ and fl S DC » ^ s ,oe 



(25) 



The behavior of Jsdc is apparently dominated by the 
smooth background term Joe near ey. As a consequence, 
the Jsdc factor in Eq. (O) no longer contributes to the 
value of A. From Eqs. (|10|), ( |TT| ) and ([l6]), one gets that 



ji = A = 2^, 



(26) 



instead of Eq. (|19|). As for the amplitudes, one finds that 

(27) 



Jo = - and J s , = — > 

^2,0 



and 



" 2^ 



Vic 2 " 



S2,0 



where we have used e « e^ near e = + . If 2i/ is not an 
integer, then the best way to evaluate the exponents fj, 
and v is, in principle, to differentiate Jsdc and Qsdc 
with respect to e and to analyze the corresponding diver- 
gences after certain orders of differentiation. The scaling 
relations fi = A = 2^, hence, provide a very strong test 
for the e — e — er assumption. If fi = 1v < 1, then the 
slope of Jsdc is negative in the range of < e < e x with 
e x = [mJ^o/J)] 1 ^ 1- ^- The observation of a negative 
slope of Jg dc near ex will apparently increase the valid- 
ity of scenario (B). But if the value of e x is very small, 
such an observation may not be practical at present. 



IV. NUMERICAL SOLUTIONS AND DATA 
ANALYSES 

We now present our numerical studies of SDC with 
the GSH equations. The numerical method for solving 
the GSH equations is based on the work by Bj0rstad 
et al. p5[ . Following Ref. 0, we choose g m = 50, 
a = 1.0 and c 2 = 2.0 for parameters in Eqs. (|l|) 
and (^). In our simulation, we take a square cell of 
size L x = L y — 1287T, which corresponds to an as- 
pect ratio r = 64. Uniform square grids with spacing 
Ax = Ay = ir/4.0 have been used, so the total num- 
ber of nodes is 512 x 512. We use the rigid boundary 
conditions i/j\b = n • VV>|b = C|s = n ' ^C\b = in 
the simulation. Here n is the unit vector normal to the 
boundary, say B, of the domain of integration. We take 
two different routes to systematically study the transi- 
tions between parallel roll states and SDC states: (A) 
We increase the control parameter e from e = 0.05 to 
e = 0.6 with steps of Ae = 0.05. We call this the roll 
branch. (B) We decrease e from e = 0.8 to e = 0.4. We 
call this the SDC branch. For e = 0.05 or 0.8, we choose 
initial conditions C( r ^ = 0) and ip(r,t = 0) as random 



variables, obeying a Gaussian distribution with a zero 
mean and a variance of 0.001. For other subsequent e's, 
we take the final results from the previous e as our ini- 
tial conditions. For each e, we wait about four horizontal 
diffusion time th before collecting data which, we hope, 
is sufficient to pass the transient regime. We run for an 
additional 20.8^,, to collect 20 instantaneous profiles for 
SDC states or 16^ to collect 10 profiles for parallel roll 
states during each data collection. 




FIG. 2. Instantaneous patterns of i))(v,t) and u> z (r,t). 
Dark regions correspond to %p(v,t) > or ui z (r,t) > while 
white regions to %l)(r,i) < or u> z (r,t) < 0. (a) tp(r,t) and 
(b) ui(v,t) at e = 0.2 on the roll branch; (c) ij)(r,t) and (d) 
u>(r, t) at e = 0.65 on the SDC branch. 

The patterns we observed are very similar to those 
found in real experiments |5||| , which can be summarized 
as the following: (A) Within the roll branch: Straight 
parallel rolls are observed at e = 0.05 up to 0.4 with a 
few defects at the boundary. Starting at e = 0.45 up to 
0.6, the rolls start to bend and focal singularities start to 
appear near the boundary. Weak time-dependence sets 
in at e = 0.45 in which the nucleation of defects along the 
sidewall occurs. Two typical shadow graph images, one 
for the order parameter tp(r,t) and another for the ver- 
tical vorticity u) z (r,t) = —W 2 (^(r,t), of the roll state at 
e = 0.2 are shown in Fig. 2. (B) Within the SDC branch: 
SDC states are observed at e = 0.8 to 0.55, whose be- 
havior has been described in detail before |qj^,^q|. At 
e = 0.5, only a few spirals still exist which mix with a 
background of locally curved rolls. Finally at e = 0.4, 
the pattern looks much like a roll state with a few de- 
fects and dislocations. In Fig. 2, we plot two typical 
shadow graph images, again one for ip(r,t) and another 
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for oj z (r,t), of SDC at e = 0.65. While the vertical vor- 
ticity u> z (r) in a roll state is almost zero everywhere, the 
corresponding field has a much richer structure in SDC. 
This suggests choosing uj z (r,t) to be an order parame- 
ter in distinguishing a roll state from a SDC state Jl5| . 
It is interesting to point out that patterns in the inter- 
val e = 0.4 — 0.6 depend on their earlier histories, i.e., 
whether they are on the roll branch or the SDC branch; 
see Fig. 3. So a hysteretic loop exists when one follows 
the two different routes. This is consistent with exper- 
imental observations that different patterns evolve from 
different initial conditions at the same e @. Owing to 
the existence of hysteresis, the transition temperature ct 
between parallel roll states and SDC states cannot be de- 
termined precisely in our study. Our rough estimate is 
e T ~ 0.45. 



(a)"" " " ' (b) 
FIG. 3. Instantaneous patterns of ip(r,t) at e = 0.55: (a) 
on the roll branch and (b) on the SDC branch. Dark regions 
correspond to ip(r,t) > while white regions to ip(r,i) < 0. 

In order to determine the character of the transition 
between rolls and SDC, we have calculated, from our nu- 
merical studies, the structure factor S(k), the two-point 
correlation length £2, the convective current J, the vor- 
ticity current SI, and the spectra entropy 5 for both roll 
states and SDC states. The numerical methods used in 
calculating these quantities are the same as in Ref. |L5) . 
The numerical uncertainties are taken as the variances 
of our data. Considering that we can, at most, take just 
a few samples of the strongly fluctuating instantaneous 
quantities, presumingly obeying Gaussian distributions 
near their corresponding time-averaged values, we believe 
that the probabilities, and hence the uncertainties, for us 
to obtain the truly time-averaged values are determined 
by those variances. From the data for S(k) of the SDC 
states, we have verified the existence of the sca ling form 
@ within our numerical uncertainties for SDC fll4[| . The 
results for £2, J, fi and S are plotted in Figs. 4-7. The 
most striking feature in these figures is the hysteretic 
loops formed by the two branches. It is also noticeable 
that the uncertainties on the SDC branch are generally 
larger than those on the roll branch, which, presumably, 
is due to the chaotic character of SDC. 

We fit the data on the roll branch with power-law be- 
haviors. To allow for the possibility that roll states might 



be unstable or metastable for e > 0.45, only data for 
0.05 < e < 0.45 are used in actual fittings, (i) We first 
use the nonlinear \ 2 method to fit the convective current 
with J = J (e — e c )^ and find that J = 0.6554 ± 0.0002, 
(j, = 1 .0054±0.0004 and e c = 0.002. The actual accuracy 
in our results may not as good as indicated. The fitting 
error for e c is very small. This e c is the measured on- 
set from conduction to convection, whose positive value 
is most likely due to finite-size effects pq . Apparently 
J on the roll branch has a mean-field exponent. The 
amplitude Jo is also in good agreement with Eq. 
provided that (x 2 ) x > 1 is not too big. (ii) Using the x 2 
method, we fit the data for the correlation length with 
£2 = ^2,o(e — e c)~", which leads to £2,0 = 13-0 ± 0.2 and 
v = 0.54±0.01. So £2 on the roll branch also has a mean- 
field exponent, (iii) Using the \ 2 method, we fit the data 
for the vorticity current with SI = Slo(e — e c ) A and find 
that Sl = (1.952 ± 0.008) x I0~ 9 and A = 2.461 ± 0.003. 
Again, the actual accuracy in our results may not as good 
as indicated. This behavior of SI is not easy to under- 
stand. The amplitude equations coupled with mean-flow 
predict, for rigid-rigid boundaries, that SI ~ e 7 / 2 for al- 
most perfect parallel rolls and SI ~ e 3 for general pat- 
terns None of these can explain the behavior of SI 
we found, which seems to be consistent with SI ~ e 5 / 2 . 
(iv) The analysis of the spectra entropy is most difficult 
since there is no theory whatsoever to describe its behav- 
ior. We simply fit it to a form H = S/, + So(e — e c ) s for the 
roll states. We find from the nonlinear \ 2 method that 
S h = 4.24±0.04, S = 2.5±0.1 and S = 1.00±0.08. The 
background term is much larger than the correspond- 
ing value for perfect parallel rolls 3& = In 2 = 0.6931. 
This, presumably, is due to finite-size effects and/or lim- 
ited computing time. The original data and their cor- 
responding fitting curves for these global quantities are 
plotted in Figs. 4-7. 

The analyses of the data on the SDC branch must 
be treated with caution given that there are three pos- 
sible scaling scenarios as discussed in Sec. III. We 



first fit the data of £2 and SI to power laws such as (i) 
£2 = 62,0 (e - e y and (ii) SI = fi Q (e - e )\ where 
60 = e c [the onset temperature in a finite system] in 
scenario (A) or eo = ex [the roll-to-SDC transition tem- 
perature] in scenario (B). Then we fit the data of J in 
accordance with our theoretical result in Eq. (|l0|), i.e., 
(iii) J = Jo(e — e c ) — ^£2" 2 > where we take the corre- 
sponding fitting results for £2 in (i). Clearly the value of 
eo is essential to determine which of the three scenarios 
is valid. In the following, we apply three different fittings 
for three values of eo and check the consistence of our nu- 
merical data against the theoretical results in Eqs. ( |To| ) 
and®. 

(a) To check whether scenario (A), i.e., e = e, is valid, 
we fix eo = e c = 0.002, whose value is given by the fitting 
of J on the roll branch. Since it is probable that SDC 
states are unstable or metastable at e — 0.4 and 0.45, 
the corresponding data might deviate from their "real" 
values in order to form the hysteretic loops. For this rea- 
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son, we disregard these two points and use only those 
data within 0.5 < e < 0.8 for fitting. We use the x 2 
method to fit our data, which leads to (i) £o = 6.8 ± 0.2 
and v = 0.72 ± 0.05; (ii) fl = (3.0 ± 0.2) x 10" 8 and 
A = 3.0±0.1; and (iii) J = 0.64±0.02and J t = 2.9±0.9. 
The original data of £2, J and fl and their corresponding 
fitting curves are plotted in Fig. 4. Apparently those 
curves fit the original data well. So scenario (A) is con- 
sistent with our numerical data. To further distinguish 
scenario (Al) \v = i] or (A2) \v > |], we find that, on 
the one hand, v — 0.72 ±0.05 from the direct fitting in (i) 
but, on the other hand, v = 0.50±0.05 from A = 3.0±0.1 
in (ii) and A = 2 + 2v in Eq. (|l9|). This discrepancy is 
likely caused by the big numerical uncertainties in our 
data. We feel that the direct fitting is more reliable and 
scenario (A2) is more likely to be true. But we cannot 
definitely rule out scenario (Al). More accurate data are 
needed to resolve this issue. 

(b) We check whether scenario (B), i.e., e = e — ct, is 
consistent with our numerical data, where the value of 
£t is determined by the fitting of £2- In contrast to case 
(a), there is no obvious reason to disregard any point 
in this scenario. So all the data within 0.4 < e < 0.8 
are used in our fittings. We first use the nonlinear \ 2 
method to fit the data of £2, which gives (i) £0 = 5.9±0.2, 
v = 0.46 ± 0.06 and e = e T = 0.27 ± 0.03. Then, we 
fix cq = 0.27 and use the \ 2 method to fit the data of 
fl and J. We find that (ii) fl Q = (8.1 ± 0.2) x 10~ 8 
and A = 2.33 ± 0.01; and (iii) J = 0.675 ± 0.001 and 

= 4.9 ± 0.1. These results are very sensitive to the 
points at e = 0.4 and 0.45. The original data of J 
and fl and their corresponding fitting curves are plotted 
in Fig. 5. The fitting of fl obviously is not good. The 
results v = 0.46 ± 0.06 in (i) and A = 2.33 ± 0.01 in (ii) 
do not satisfy the scaling relation A = 2v in Eq. (p6|). So 
scenario (B) with er — 0.27 is unlikely to be true. 

(c) We check whether scenario (B), with et determined 
by the fitting of fl, is consistent with our numerical data. 
As in case (b), we use all the data within 0.4 < e < 0.8 in 
our fittings. We first use the nonlinear x 2 method to fit 
the data of fl and find that (ii) fl = (4.5 ± 0.3) x 10~ 8 , 
A = 1.41 ± 0.07 and e T = 0.348 ± 0.005. Then we fix 
eq = 0.348 and apply the \ 2 method to fit the data 
of £2 and J. We find that (i) £2,0 = 6.55 ± 0.09 and 
v = 0.295 ± 0.008; and (iii) J = 0.667 ± 0.001 and 

= 4.6 ±0.1. Again, these fitting results are very 
sensitive to the points at e = 0.4 and 0.45. The origi- 
nal data of £2, J and fl and their corresponding fitting 
curves are plotted in Fig. 6. Notice that, since v in (i) 
is less than ~, the slope of J is negative in the range of 
< I = e - e T < e x with i x = [2vJ ( /^ fi J ] 1/(1 ~ 2l/) . 
In the present case, one has e x = 3.20 x 10~ 3 , which 
perhaps is too small to be checked by real experiments 
or simulations. From the pure-data-fitting point of view, 
the fittings in this case are as good as the fittings in case 
(a). So, without the benefit of our theoretical results, 
scenario (B) with et = 0.348 could be (incorrectly) ac- 



cepted. But, the exponents v — 0.295 ± 0.008 in (i) and 
A = 1.41 ± 0.07 in (ii) are not even close to satisfying the 
scaling relation A = 2v in Eq. (|2^). So this scenario can 
be ruled out by our theory. From this, together with the 
discussions in case (a) and (b), we conclude that scenario 
(B) is unlikely to be valid for SDC and scenario (A) is 
consistent with our numerical data and our theory. 




e 

FIG. 4. Plots of £^ 2 , J and fl vs. e for both parallel rolls 
and SDC. The corresponding fitting curves for SDC are de- 
scribed in (a) in Sec. IV. The labels are used as the follow- 
ing: SDC — the numerical data with error bars on the SDC 
branch; (a) — the fitting curves for SDC described in (a) in 
the text; Roll — the numerical data with error bars on the 
roll branch; (p) — the fitting curves for parallel rolls. To indi- 
cate the parallel roll states for e > et and the SDC states for 
e < £t may be metastable, the corresponding fitting curves 
in (a) and (p) are plotted with dotted lines. The transition 
temperature is estimated roughly at er = 0.45. 

It is worthwhile to mention that, in all cases above, 
the value of Jo agrees with the theoretical result 2/3 in 
Eq. ([l(]); the value of is also consistent with the the- 
oretical prediction > 8/3. On the other hand, the 
theoretical results in Eqs. (|2l|) and (ps|) predict that 
fl = 7.4 x 10~ 4 in case (a), ft = 7.2 x 10~ 5 in case (b), 
and flo = 9.7 x 10~ 5 in case (c). All these predictions 
are several orders larger than the corresponding numeri- 
cal results. The reason for such big discrepancies is not 
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clear to us. 

We have also studied the behavior of the spectra en- 
tropy. Owing to the lack of any theory, we simply fit 
the data of 3 to a form 3 = 3& + So(e — eo) A for the 
SDC branch. We apply two different fittings: (a) We 
fix eo = e c and 3& = 4.24 [from the fitting for the roll 
branch], and use the x 2 method to fit those data within 
0.5 < e < 0.8, which leads to 3 = 4.79 ± 0.03 and 
5 = 0.18 ± 0.02. (b) We fix 5 fc = 4.24 and use the 
nonlinear \ 2 method for all data within 0.4 < e < 0.8, 
which gives that 5 = 5.13 ± 0.05, 5 = 0.12 ± 0.01 and 
eo = 0.37 ±0.02. We have also tried other alternative fit- 
tings such as using the nonlinear y 2 method for all data 
within 0.4 < e < 0.8 and fixing e = 0.27 [from (i) in (b)] 
or e = 0.348 [from (ii) in (c)], but none of them gives a 
reasonable fit. The fitting curves in (a) and (b) and the 
original data of 3 are plotted in Fig. 7. At this stage, 
the behavior of 3 is the most unclear one among all the 
time-averaged global quantities defined in Sec. O. 




FIG. 5. Plots of £^~ 2 , J and Q. vs. e for both parallel rolls 
and SDC. The corresponding fitting curves for SDC are de- 
scribed in (b) in Sec. IV. The labels are used as the follow- 
ing: SDC — the numerical data with error bars on the SDC 
branch; (b) — the fitting curves for SDC described in (b) in 
the text; Roll — the numerical data with error bars on the 
roll branch; (p) — the fitting curves for parallel rolls. 




e 

FIG. 6. Plots of £^ 2 , J and Q vs. e for both parallel rolls 
and SDC. The corresponding fitting curves for SDC are de- 
scribed in (c) in Sec. IV. The labels are used as the follow- 
ing: SDC — the numerical data with error bars on the SDC 
branch; (c) — the fitting curves for SDC described in (c) in 
the text; Roll — the numerical data with error bars on the 
roll branch; (p) — the fitting curves for parallel rolls. 



V. DISCUSSION AND CONCLUSION 

In the previous section, we conclude that scenario (A) 
is valid for SDC. This means that £2 diverges at e = 0, 
and J and vanish at e = 0. At first sight it seems 
puzzling that all properties of SDC are controlled by e, 
instead of e — ex- To understand this, we propose an ex- 
planation for this scenario, which is somewhat similar to 
that in the hexagon-to-roll transition in non-Boussinesq 
fluids (20). In the latter case, the transition from hexag- 
onal states to parallel roll states occurs at finite ex- Al- 
though the roll attractor is unstable for small enough e 
and metastable against the hexagonal attractor for even 
slightly larger e, the properties of the parallel roll states 
are all controlled by e, not by e — ex [^0|. Clearly one can 
imagine a similar picture for the roll-to-SDC transition. 
While the SDC attractor seems to be either unstable or 
metastable against the roll attractor for sufficiently small 
e, as an intrinsic convective state, the properties of SDC 
are controlled at the conduction to convection threshold, 
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not where it starts to emerge as the stable state. The 
existence of two different attractors has been suggested 
by experiments ]lT||I^ ]. The basins and the stability of 
these two attractors are still unclear at present. 
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FIG. 7. Plot of S vs. e for both parallel rolls and SDC. 
The different fitting forms (a) and (b) for SDC are discussed 
in the text. The labels are used as the following: Roll — the 
numerical data with error bars on the roll branch; SDC — 
the numerical data with error bars on the SDC branch; (a) 

— the fitting curve for SDC described in (a) in the text; (b) 

— the fitting curve for SDC described in (b) in the text; (p) 

— the fitting curve for parallel rolls. To indicate the parallel 
roll states for e > er and the SDC states for e < er may 
be metastable, the corresponding fitting curves in (a) and (p) 
are plotted with dotted lines. The transition temperature is 
estimated roughly at er = 0.45. 

The establishment of scenario (A) indicates that the 
transition between the parallel roll states and the SDC 
states is first-order. This conclusion is also supported by 
the following: 

(1) Our theory predicts discontinuities in the value and 
the slope of J at [see the discussion following Eq. 
(|l5|)]. This is a typical signature of a first-order transi- 
tion. 

(2) The presence of hysteretic loops in Figs. 4 - 7 is 
a strong indication of a first-order transition. A differ- 
ent hysteretic loop has also been reported by others for 
the GSH model [|13|. Although it is arguable that hys- 
teretic loops may be found in a second-order transition 
if the computing time is not long enough to overcome 
the effects of critical slowing-down (which occurs when 
the correlation time t approaches infinity), it is doubtful 
that the loops in that case can be as distinctive as what 
we found in Figs. 4-7. 

(3) As we described in Sec. |v|, the convective patterns 
depend on the processes leading to them, which is also 
observed in experiments |Q . This fact suggests that the 
two competing attractors are either both stable or one is 
stable while the other is metastable for some positive e. 
Such a stability property is typical in first-order transi- 
tions. On the contrary, in second-order transitions one 



of the two attractors should change from stable to unsta- 
ble while the other changes from unstable to stable as e 
moves across er- 

(4) From Figs. 5 - 7, it is easy to see that, if sce- 
nario (B) is valid, then the fitting curves of all the time- 
averaged global quantities on the SDC branch will cross 
those on the roll branch. But there is no evidence from 
our numerical calculation supporting such a crossing. 

As we mentioned in Sec. [|, an earlier experiment with 
a circular cell Q found that £ 2 diverges at e = with 
a mean-field exponent, while the correlation time r ei- 
ther diverges at e = with a non-mean-field exponent 
or diverges at the roll-to-SDC transition temperature ex 
with a mean-field exponent. However, a recent experi- 
ment with a square cell (lj] concluded that £2 diverges 
at ex with a very small exponent v. We now comment 
on these experiments and our study. 

Regarding our numerical study, we cannot rule out 
that the roll-to-SDC transition in the GSH model has 
a different character from those in real experiments al- 
though this seems unlikely. We also cannot rule out that 
our numerical solutions are still in the transient regime 
even though we have waited for about four horizontal dif- 
fusion times th before collecting data over an interval of 
several th for each e. Furthermore, as we discussed in Sec. 
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our numerical data are not accurate enough to deter- 
mine independently which of the three scenarios is true 
for SDC. As a result, we have to rely on our theoretical 
predictions to resolve this issue. Finally we disregarded 
two data points in our analysis for scenario (A) on the 
basis that these two points deviate from their "real" val- 
ues to form hysteretic loops. This introduces a certain 
arbitrariness in determining which data points deviate. 
These shortcomings in our numerical study weaken the 
validity of our conclusion. 

We notice that the data of the correlation length £2 
for parallel rolls and SDC were analyzed together in the 
experiment by Morris et al. Jl~T[ ] . which we think is not 
justified. Considering that the parallel roll states and 
the SDC states are intrinsically different, we believe it is 
necessary to separate their data in the analysis, such as 
we did in Sec. |y|. Such a separation was implicit for the 
data of the correlation time t since t = +00 for steady 
states such as parallel rolls and, in principle, only the 
data for SDC are available. A divergence at e = with a 
non-mean-field exponent was found to be consistent with 
the data of r for SDC O], It is not clear to us whether 
a similar conclusion can be reached for £2 if its data for 
SDC are analyzed separately. 

The conclusion by Cakmur et al. jjJI that £2 diverges 
at ex with a small exponent v is different from that in 
Ref. Q and ours. These authors suggested that the dif- 
ferent geometry of experimental cells is the cause for the 
different results between theirs and that in Ref. Jll|. If 
so, then it is not clear to us what the real behavior is 
in an infinite system. In our calculation, we also used 
a square cell. We think that there may exist a different 
interpretation of the experimental data for £2 obtained 







in Ref. p2[. We believe that, to convincingly show that 
£2 for SDC diverges at ex, one must have sufficient data 
points whose £2 are much larger than the corresponding 
typical values of parallel rolls (for e away from both 
and ex)- Otherwise, the data for £2 for SDC may simply 
approach those for parallel rolls to form a hysteretic loop 
near ex, instead of diverging at ex- Since no data for £2 
for the parallel roll states was given and the number of 
data for SDC near ex was relatively limited, we feel the 
conclusion by Cakmur et al. |l2j| should be treated with 
caution. 

As we discussed in Sec. IV, our theory plays an im- 
portant role in determining the nature of the roll-to-SDC 
transition. So it is very important to check our theo- 
retical predictions by real experiments. One important 
prediction by our theory is that there exist discontinu- 
ities in the value and the slope of J at ex- However, 
we realize that no discontinuity in J has been reported 
by experiments. The reason for this is not clear to us. 
We conjecture that finite size effects might play a role. 
From Fig. 1, we find that the discontinuity is larger for 
smaller Prandtl number a. So it would be interesting to 
see whether experiments can confirm or rule out such a 
discontinuity in J by using a small a. Another impor- 
tant prediction from our analysis is the behavior of the 
time-averaged vorticity current f2. Since direct measure- 
ments of seem to be very difficult in real experiments 

@, we think it valuable to calculate ft by solving Eq. 
[the corresponding version before rescalings can be 
found in Ref. |14|] or its improved versions numerically, 
with the experimental results of ip(r,t) as input. Such 
a calculation will not only help to clarify the nature of 
the roll-to-SDC transition, but also provide an additional 
experimental test on our theory jlj]. It would also be 
useful to calculate the time-averaged spectra entropy as 
suggested in Refs. [jl2|,[l5), even though there is no theory 
to predict the behavior of this quantity. 

In summary, we conclude from our numerical studies 
and our theoretical results that the roll-to-SDC transition 
is first-order in character. We found that the correlation 
length £2 for SDC diverges at e = 0, not at the transition 
temperature ex- However, since the uncertainties in our 
data are unpleasantly large and the data points we have 
are unsatisfactorily few in number, we cannot determine 
definitely whether the exponent of £2 is mean-field or not. 
So further investigations are necessary to draw a definite 
conclusion. In this regard, a theoretical calculation of £2 
for SDC is highly desirable. A theory to describe the roll- 
to-SDC transition is essential. Finite size effects should 
also be studied carefully. 
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